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Abstract 

Models relating to the Species-Area curve usually assume the existence of species, 
and are concerned mainly with ecological timescales. We examine an individual- 
based model of co-evolution on a spatial lattice based on the Tangled Nature model 
in which species are emergent structures, and show that reproduction, mutation and 
dispersion by diffusion, with interaction via genotype space, produces power-law 
Species-Area Relations as observed in ecological measurements at medium scales. 
We find that long-lasting co-evolutionary habitats form, allowing high diversity 
levels in a spatially homogenous system. 
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i 1 Introduction 



2 The number of species in a given region can be seen as a product of the 

3 evolutionary history of speciation, extinction and migration to that region. 
Time variations in an ecology, whether induced by population dynamics or 

5 evolutionary dynamics, are caused by processes operating at the level of indi- 
e viduals; taxonomic structures, like species and genera, are emergent entities 
7 produced by the unceasing action of reproduction, mutation and annihilation 
s of individuals. Hence it should be possible to derive the stability properties, 
s abundance and distribution of species from a 'microscopic' description in terms 
of dynamics at the level of individual organisms. Such a framework must be 

11 able to act as a unified explanation of ecological structures such as the Species 

12 Area Relation (SAR) and the Species Abundance Distribution (SAD) together 
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13 with evolutionary aspects such as the temporal variation of the macroscopic 

14 averaged extinction rate and intermittency in the extinction events. 

15 The relationship between the number of species observed in an area and the 

16 area's size is one of the most basic questions in ecology but it is still the 

17 subject of much debate. The number of species found in an area could increase 
is with area size simply because more individuals are counted, and the form of 

19 this relation may be very different depending on the counting method used 

20 and details of the area [1][2]. For most measurement scales on non-island 

21 systems it seems that a power-law - (diversity) oc (area) 2 - may be an accurate 

22 description, for the majority of fauna and flora types. However, for other scales 

23 and for some data, other forms have been successfully fit [3]. Here we consider 

24 those systems for which a power law provides a good fit - we will comment 

25 below on possible effects not included in our model which may be responsible 

26 for observed non-power law forms. 

27 Dynamical models typically assume the existence of a set of species as given 

28 structures classifying individuals. The dynamics at the individual level then 

29 determines how the assumed species are, say, populated and spatially dis- 

30 tributed. Particularly impressive examples of this type of models are Hubbell's 

31 [4] 2001 'Unified Neutral Theory' and Durrett and Levin's [5] 1996 spatial voter 

32 model. In the neutral models [4] [5] [6] [7] all individuals have the same birth, 

33 death and migration rate independent of which species they belong to. Sole, 

34 Alonso and McKane [8] have introduced a more general set of models in which 

35 an interaction matrix allow the assumed set of species to vary in their proper- 

36 ties. Choosing specific forms for the interaction matrix reduces this model to 

37 a number of previously considered models - among these is Hubbell's neutral 

38 model. Realistic SAD and SAR are obtained from these models even in the 

39 case of neutrality between species. The SAD and the SAR has in addition 

40 been explained by an attractive geometric approach by Harte and co-workers 

41 [9], who replaced dynamics by the assumption that the spatial distribution of 

42 species is self-similar and fractal; a prediction which has been confirmed from 

43 field data on birds in the Czech Republic [10]. They concluded that a power 

44 law SAR was equivalent to a community level fractal distribution of species. 

45 The Tangled Nature model (TaNa for short) is an attempt of developing a log- 

46 ically simple approach to evolutionary ecology. From a few fundamental and 

47 generally accepted microscopic assumptions, macroscopic phenomenon such 

48 as macroevolution and ecological structures emerge. The model is individual 

49 based with fluctuating population size, and the mutation prone reproduction 

50 occurs with probabilities determined by the interaction between co-evolving 

51 organisms. The long time macroevolution in the model is consistent with ob- 

52 served temporal characteristics [11], the SAD compares well with observation 

53 [12] and most recently the model has been used to understand microbiological 

54 experiments on the relation between diversification and interaction [13] . In the 
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55 present paper we demonstrate how the Tangled Nature approach can be used 

56 to understand the SAR from an evolutionary individual based view point. 

57 We will be introducing spatial aspects into the non-spatial TaNa model, in 
ss order to measure the SAR. Essentially all good dispersion models produce 
59 reasonable fit with data (usually a power law) - e.g. the spatial models dis- 
eo cussed above. Power-laws are often observed in field data, but not universally 
ei [14], and we hope to eliminate two of the possible causes of the deviation - 

62 interactions and localisation (i.e. deme structure). The interaction permitted 

63 in our model provides approximate power-law SAR regardless of strength, so 

64 inhomogeneity in migration or resource is a more likely source of observed 

65 deviations from power law in real systems, as such inhomogeneities are not 
ee considered here. 

e/ Here we consider species as dynamical quantities that emerge in genotype 

es space. We allow for spatial extension in a homogeneous physical environment, 

69 breaking the population into a number of spatial locations (with each species 

70 type forming separate demes) which in our model permits the construction 

71 of co-evolutionary habitats 1 of interacting species within each lattice point. 

72 Individuals move by random dispersion as in the models mentioned in the 

73 previous paragraph. The co-evolutionary habitats survive for very long time 

74 periods, during which local species abundances fluctuate around some average 

75 level. Inside these habitats equivalence of individuals is observed, as a result of 

76 adaptation. The offspring probability of an individual depends on its genotype 

77 and on the composition of the local community in the local genotype space. 

78 All individuals are subject to the same annihilation rate and only individuals 

79 that have evolved genotypes with an offspring probability that matches the 
so killing probability are able to constitute species with a degree of temporal 

81 stability. This leads to a certain degree of equivalence or neutrality to emerge 

82 amongst the dynamically generated species. Since the offspring probability of 

83 an individual depends on the local occupation of genotype space, when indi- 

84 viduals disperse to other habitats they usually do not have the same offspring 

85 probability as the members of the habitat they enter. If species composition 
se begins to change locally, then the entire habitat is usually affected, disrupting 
87 the local species composition. 

ss Interaction allows for the extinction of well-established species on ecological 

89 timescales in the right invasion circumstances, giving realistic Species Abun- 

90 dance Curves (approximately log- normal [11]). Although species in the Tan- 

91 gled Nature model are dynamical and emergent, properties associated with 



We use the term co-evolution in the weak sense of species that have adapted due to 
interactions with other species. We will also refer to these 'co-evolutionary habitats' 
as simply 'habitats' for brevity, as they are the only kind of habitat possible in our 
model. 
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92 random dispersal such as power-law SAR are observed. The interaction allows 

93 distinct species to be separated in genotype space, in contrast with neutral 

94 models. In hypercubic genotype space and in the absence of interaction species 

95 are clustered around a mean with separation occurring only by fluctuation and 

96 persisting for short timescales [15] (this also tested for the non-interacting ver- 

97 sion of our model, where the population essentially moves stochastically as one 

98 coherent cluster through genotype space). 

99 The original Tangled Nature model defined by Christensen et al. [11] has no 

100 spatial component, which we introduce by running copies of the model con- 

101 currently on a square lattice, allowing for interaction between lattice-points by 

102 migration. The interaction between individuals at adjacent sites is therefore 

103 indirect, acting through genotype space only via the distribution of migrants, 

104 and the spatial aspect is discrete. However, we can easily compare our results 

105 to that of the original model which has stability properties known to be close 
we to observed systems [16] [17] [12]. The motivation for our approach is that gen- 
107 era that can move (animals and bacteria), or whose offspring can compete 
los over distance for space (most plants) are modelled as locally well mixed, with 
109 spatial aspects considered on larger scales. 

no We begin with a recap on the non-spatial Tangled Nature model and its major 

in features. Then we detail our simple extension to the model introducing spatial 

n2 dimensions using a square lattice of models. 



n3 2 Definition of the Model 

H4 We now define the Tangled Nature model. We will be constructing the model 

us on a periodic square lattice of length X. Specific points on the lattice are 

lie referred to by their co-ordinates (x, y) . Each point on the lattice may contain 

H7 any number of individuals who, on any given time step, may migrate with 

us probability p m0 ve to a neighbouring lattice point (our neighbourhood includes 

n9 diagonals, and therefore is 8 lattice-points). On each lattice point we run a 

120 TaNa summarised below and described in [11] [17], with interaction between 

121 lattice points via migration. Each lattice point contains a number of species, 

122 made up of explicitely modelled individuals. Similar approaches have been 

123 used many times, e.g. with each lattice point containing a local food web [18], 

124 or being used as the basic unit instead of individuals for models in which the 

125 two scales can be well separated (Gavrilets book [19] considers this and many 

126 other situations). Such separation of scales is not possible in our model, as 

127 the specifics of individuals control the invadability and stability of the local 

128 population. 

129 The Tangled Nature model represents individuals as a vector S Q = S%, SI) 
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130 in genotype space. The Sf take the values ±1, and we use L = 20 throughout. 

131 Each S a string represents an entire species with unique, uncorrelated interac- 

132 tions, i.e. genotype space is coarse-grained. The small value of L is necessary 

133 for computational reasons as all genotypes exist in potentia and have a des- 

134 ignated interaction with all other possible organisms. It is also possible to 

135 define the model slightly differently in terms of smooth traits, and correlate 

136 interactions over the trait space [20]. 

137 We refer to individuals by Greek letters a, (3, ... = 1,2, ...,N(t) for a specific 
iss lattice point (x, y). Points in genotype space are referred to as S a , S b , and 
139 many individuals (from any real-space location) may belong to a point in 
wo genotype space S a . 

141 In the TaNa model, all individuals are considered to die with equal probability 

142 Pkui, so it is most appropriate to systems where competition is for offspring 

143 space or resources (plants or bacteria, for example). Only the probability to 

144 produce offspring is controlled by their interactions; however, the model is 

145 qualitatively the same regardless of whether varying killing or reproduction 
we rates are used [11]. Reproduction occurs asexually, and on a successful repro- 
147 duction attempt a daughter organism is produced which will be mutated with 
us probability p mu t- When an individual a is chosen for processing, it will repro- 
149 duce with probability: 



PoffiS '^-l + exp^SV)] G(0 ' 1) (1) 

p ff is defined in this way as it is the simplest way to translate H(S a ,t) into 
a reproduction probability. H(S a ,t) is defined in Equation 2 and contains 
the bulk of the model, consisting of interaction and competition. It is the 
average interaction (first term) and resource competition (second term) with 
all other individuals in the same spatial location. Interactions are considered 
as an average (hence dividing by the population size N(t)) and we write it 
as a sum over all species rather than individuals, as individuals of the same 
species are identical. 



H(8°, t) = —L- • £ J(S°, S)n(S, t) - »N(t) (2) 



c is a parameter controlling the interaction strength, N(t) is the total num- 
ber of individuals at time t and n(S,t) is the number of individuals with 
genotype S at that point, /i controls the carrying capacity of the system, pre- 
venting population growth when iV is of the order \j \x. The interaction matrix 
J(S a , S) represents all possible couplings between all genotypes, each gener- 
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164 ated randomly in the range (—1,1), being non-zero with probability O. Since 

165 the functional form of J(S a , S b ) does not affect the dynamics, provided that 
lee it is non-symmetric with mean 0, we choose a form of the interaction matrix 

167 that speeds computation [11]. In the spatial version, we use the same S but 

168 allow the individuals to be located at a point in space, such that a = a(x,y), 

169 N = N(x, y, t) and n = n(x, y, S, t). 

no Since the elements of J are generated randomly, the pairwise interactions 

171 can be of the following types: mutualism (both positive), competition (both 

172 negative) and predator/prey (or parasitic) relations (one positive and one 

173 negative). We do not allow for one-way interactions such as amensalism, apart 

174 from in the case where one interaction is randomly generated to be very small. 

175 Also note that even in the case of extreme mutualism, resource is limited and 

176 competition will occur as the population increases, and so the negative term 

177 fJ,N(t) in Equation 2 is large. 

178 The interactions modelled here are very general, though must occur through 

179 some medium which is not modelled explicitly. For bacterial systems, this 

180 would be in the form of chemicals, meaning that the resource is modelled 

181 to some degree, but for plants it is more likely to be direct competition for 

182 offspring space. The limiting resource, controlled by \x is different to any inter- 

183 action facilitating resource, and might be space or a food source depending on 

184 the system under comparison. There is only one 'type' of resource, however, 

185 and as such we are only really modelling within a single trophic level, amongst 
we individuals concerned with the same basic resource. Thus, our model can be 
187 compared with data for herbivorous birds, or bacteria, or crop plants, but 
las only for a predator-prey system when individuals on different trophic levels 

189 still compete for space. This is not a problem for this papers purposes as most 

190 SAR data is drawn from a single family of species. We are trying to model 

191 both the obvious food-web interactions as well as the multitude of perhaps 

192 weaker, hidden interactions. It would be simple to add a number of additional 

193 resource types, with species drawing variously from different resources, but 

194 this adds a level of complexity unnecessary for the current questions. This is 

195 instead considered as an extension to the model [20]. 

we In an offspring individual, each Sf is mutated (flipped from 1 to -1, or from -1 

197 to 1) with probability p mut from the parental Sf . Thus mutations are equiva- 

198 lent to moving to an adjacent corner of the L-dimensional hypercube in geno- 

199 type space, as discussed in [11]. 

200 A time-step consists of choosing a spatial lattice point with probability pro- 

201 portional to the population of the lattice point N(x, y, t). Then an individual 2 

2 In previous versions a different individual was chosen for reproduction and killing 
actions. Here we select only one individual and process it for reproduction, killing 
and movement for code efficiency reasons - above the level of fluctuations the two 
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202 a is chosen randomly from that lattice point. 

203 • a is allowed to reproduce with probability p ff- 

204 • a is killed with probability Pkui- 

205 • If the killing attempt was unsuccessful, a is moved to an adjacent lattice 

206 point with probability p moV e- Thus the effective p e JJ me = (1 - Pkiii)Pmove- 

207 We define a generation as the amount of time for all individuals to have been 

208 killed, on average, once. For a stable population size, this is also the time 

209 for all individuals to have reproduced once, on average. Generations therefore 

210 are overlapping, and individuals have an exponential lifetime. The choice of 

211 constant pkui does not appear to affect the general results - if we reversed 

212 the situation and allowed constant p Q ff whilst varying fitness via Pkui, the 

213 same behaviour is observed (as the equilibria has p a ff ~ pkm for all species). 

214 We should therefore not observe results that are specific to either high infant 

215 mortality or high adult competition mortality, but we should observe features 

216 common to both competition types. 

217 Although our model is asexual, we are operating on a sufficiently course- 

218 grained level that sexual reproduction can be considered as only possible be- 

219 tween two individuals of the same genotype, and therefore is identical to the 

220 asexual case in our model, apart from when the abundance of a species is so 

221 low it would not be able to find a mate. Whilst this permits comparison with 

222 data from both sexually and asexually reproducing species, this approxima- 

223 tion will be invalid for many cases; we do not consider cross-over effects, for 

224 example. One can think of our genotype space as modelling the genes that 

225 effect fitness, with 'neutral' variation permitted in a type without being ex- 

226 plicity modelled. Some of the effects of sex could be incorporated into the 

227 mutation probability - others must simply be ignored. We have not yet found 

228 any population level data that significantly contradicts our model, although 

229 clearly we miss a lot of the fine detail. A discrepancy between our model and 

230 observed data which is only present for sexual species could shed light on the 

231 population level effects of sexual reproduction. 



232 3 Behaviour of the model 

233 We will first review the behaviour of an isolated system, and then use this to 

234 help interpret the results on an X by X square lattice with periodic boundary 

235 conditions. 

236 Unless otherwise stated, the parameters used will be: 6 = 0.25, c = 0.05, \i = 
methods are equivalent. 
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237 0.05 j Vmut — 0.01 and pkui — 0.1; see [11] for more details. These are chosen to 

238 keep the population of the entire system from exceeding about 30000, keeping 

239 computation to reasonable levels and allowing for averaging. The population 

240 of a specific lattice point is low compared with previous studies (around 300 

241 in this study), increasing the strength of stochastic effects - hence the other 

242 parameters are chosen to cancel out this effect to some degree. Although the 

243 mutation rate is unrealistically high, it still reproduces the correct qualitative 

244 effects found in real systems [17], and simply gives a higher turnover of quasi- 

245 evolutionary stable strategies as defined in Section 3.1. It should be stressed 

246 that the qualitative behaviour observed here is seen at mutation rates down 

247 to 10~ 8 . 



248 3.1 The isolated TaNa model 

249 We briefly review the behaviour of a single TaNa model as given by [11] [16]. 

250 The model exhibits a number of quasi-evolutionary stable strategies (q-ESSs) 

251 in which the frequency distribution in genotype space remains constant (with 

252 some small fluctuations); these q-ESSs are also observed in differential equa- 

253 tion style models [21] [22]. The q-ESSs are named after the Evolutionary Stable 

254 Strategies [23], or ESSs, found in game theory. If we think of competing in- 

255 dividuals which may adopt a strategy for survival, then an ESS is a strategy 

256 which, if adopted by the entire population, will not be invadable by any other 

257 strategy. The strategy of an individual defines its actions in all circumstances; 

258 in our model the strategy is the list of interactions with all other types. It 

259 is the strategy of the population as a whole that is important here, given by 

260 the proportion of individuals having each individual strategy. A stable strat- 

261 egy is thus a set of individuals who cannot be invaded by an increase in any 

262 of the other types (that is, if one type gains population, it loses interactions 

263 and therefore will lose population). However, because we include mutation, 

264 the strategy must also be stable to an influx of all local types. This list of 

265 local types is only a subset of all possible invader strategies, and so a popu- 

266 lation may be quasi-stable; that is, stable to all local mutations but not to 

267 distant genotypes which can only be reached by stochastic fluctuations (as 

268 they are separated from the population by a fitness minima). These distant 

269 genotypes can do well in the q-ESS population, and therefore destabilise it as 

270 their population grows. 

271 We operate with parameters that give a reasonable number of q-ESS switches 

272 within the first 50000 generations - most of the work analysing this region 

273 was done in [17]. During these q-ESSs (shown in Figure 1 (b)), a number of 

274 genotypes (the 'wildtypes') are highly occupied - other genotypes are only 

275 present due to mutation from the wildtypes, and are frequently eliminated by 

276 stochastic events (see Figure 1 (a)). As our genotype space is coarse-grained, 



8 



these 'sub-species' do not inherit interaction properties from a wildtype - de- 
spite this, a natural species-concept emerges as a simple result of interaction 
in a genotype space. Thus our diversity measure is the 'wildtype diversity': 
simply the number of wildtypes in the system. Wildtypes are defined as geno- 
types with occupancy of eight or greater (a definition which is valid only for 
these parameter ranges). We have tested other diversity measures such as the 
Shannon- Wiener Index and our results are qualitatively the same regardless 
of measure used, but these are primarily designed to avoid sampling problems 
[24] and so are less relevant to computer simulations. 
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Fig. 1. (a) Species abundance, or relative occupancy of points in genotype space, averaged over all 50000 
generations and 380 runs. There is a distinct difference between those genotypes with occupancy less than 
about 6 and those with an occupancy greater than 16, with only a very small amount in between, most 
of which come from transition period species, (b) An example occupancy plot showing all species with 
occupancy n a > 8 at time i as a dot. Species are not meaningfully ordered. q-ESS periods are shown as 
horizontal lines, with most transition periods (apart from the very slow one from around time 11000 to 
20000) too short to see on this scale. 

286 In [16], it is shown that the average q-ESS length increases with time, due to 

287 increasing stability in the network of active interactions, increased population 

288 size and hence increased diversity (as larger populations are more likely to be 

289 stable to stochastic fluctuations, and q-ESS interactions tend to be positive). 

290 Note that these effects occur only on average - it is possible for the system to 

291 move to a less stable, smaller population after a disordered phase, and it is 

292 also not always true that higher populations are more stable (or more diverse), 

293 just that they are on average. 

294 During the q-ESS, wildtype occupation fluctuates around some constant level, 

295 and sub-species appear and dissappear by mutation, without affecting the 

296 stability of the q-ESS state. Biologically, a q-ESS has all species in a q-ESS 

297 occupying a fitness maxima (that is, all mutations have lower fitness - fitness 

298 meaning offspring probability in this case), which the system has found dur- 

299 ing a transition. Each species in the q-ESS must have reached a population 

300 equilibrium, so that p Q ff ~ ptui, and all mutants from each species must have 

301 p ff < pkui when their own population is low. This is easier to achieve for a 

302 low diversity, but when a stable state is found at higher diversities, the chance 
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303 that an invader will destabilise the q-ESS is lower as invaders will be at signif- 

304 icantly lower fitness on average (due to the increase in the average population 

305 iV from those positive interactions). It is therefore of interest to analyse the 

306 transition more closely, in order to understand why the q-ESS forms in the 

307 way it does. 

308 Transitions appear in many forms, depending on the configuration of the geno- 

309 type space surrounding the wildtypes. There are two events that can force a 

310 q-ESS to end: 

3n • If a genotype with p Q ff > pkm can be reached, then there will be a period 

312 where the mutant population is still vulnerable to accidental extinction, 

313 followed by an exponential growth period if the mutant population grows 

314 large enough. This will usually quickly upset the configuration of the local 

315 population, leading to transition. 

316 • If one of the wildtype species had low average population then it can become 

317 accidentally extinct. In some cases other species will not depend on this 

318 species and the system enters a similar q-ESS with reduced diversity; in 

319 other cases, the stability of the q-ESS is upset and a transition occurs. 

320 Once the system enters a transition, one of the following may happen: 

321 • The disruption is minor and the system remains stable with a new q-ESS 

322 configuration. The transition period is not well defined in this case. 

323 • Wildtype species no longer all have p Q ff = Pkm- The populations will change 

324 in order to regain this relation. It is possible that a species may become 

325 extinct, leading to stage 2 above. 

326 • One of the low population mutant species in the system will gain p a ff > pkm 

327 and so will enter phase 1 above. 

328 Clearly, this is an iterative process and can last for a very long time - forever if 

329 c or Pmut are very large, so pushing the system past the 'error threshold' [17]. 

330 It is additionally complicated because these processes are all really running 

331 simultaneously, and responding to each other. What is clear, though, is that 

332 there is always favoured species in the system, and from simulations we see 

333 that the number of favoured species does not change significantly from q- 

334 ESS periods. In [11] it is shown that transition periods retain the distinction 

335 between (short lived in this case) wildtypes and mutants, resulting in a very 

336 similar (possibly identical) SAD. Since the transition periods are very short, 

337 any deviation from the q-ESS SAD is negligible and for an instantaneous 

338 observation they are indistinguishable (as stochastic noise is high). Transitions 

339 also provide a way for a species to mutate to a distantly related genotype 

340 quickly. Because there is a high interaction between all types, and the number 

341 of types is often quite high, most configurations are not q-ESS. It is therefore 

342 unlikely that the initial invaders of a q-ESS will be successful in the long 
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343 term - they instead will be in turn invaded by a second set of mutants. This 

344 process continues until a q-ESS is found, and so there is an effective selection 

345 gradient away from the wildtypes during this time, leading to very large and 

346 fast changes in genotype acting for short periods of time. 

347 The species abundances are of log-normal form as observed in many real sys- 

348 terns [12] provided that the interaction probability G is high, as in the cases 

349 we consider, and the lifetime distribution for species is wide-tailed as in real 

350 data [11] (following a power-law). More details on the network properties of 

351 the Tangled Nature model is available from [12], and an in-depth analysis 

352 of the time dependence of many of the observables such as diversity and to- 

353 tal population is presented in [16]. Similar work by Zia and Rikvold [25] [26] 

354 deals with a simplification of the non-spatial case. In both models the q-ESS 

355 wildtypes are characterised as different to transition period wildtypes because 

356 their mutants do not interact favourably with the q-ESS population, and so 

357 are suppressed. 

358 3.2 The Tangled Nature Model on a spatial lattice 

359 We now introduce a square spatial grid of length X, each containing a TaNa 

360 model, and allow the lattice-points to interact by migration; migration proba- 

361 bility refers to the chance of moving to any neighbouring site, chosen randomly 

362 from the 8 nearest neighbours, and we assume a periodic boundary. Just this 

363 simple addition to the basic TaNa model gives rise to naturally occurring 

364 Species-Area Relations, or SARs. 

365 Unlike the non-spatial version of the model, initial conditions are relevant. 

366 All possible starting configurations reduce to one of the following two initial 

367 conditions: 

368 (1) Individuals are generated with a random genotype and placed on a ran- 

369 dom lattice point until the total starting population is reached. 

370 (2) A single lattice point is allowed to evolve as a separate system until a 

371 q-ESS is formed. This q-ESS is copied to all other lattice points to give 

372 a quasi-stable, identical initial starting condition at all points. 

373 Procedure 2 represents the biological case where a small species set is exposed 

374 to a larger spatial range, and so colonises it. The initial q-ESS used in pro- 

375 cedure 2 has stability properties that can differ greatly - see Figure 2. It can 

376 vary in absolute stability (how long it will last for), but spatial duplication 

377 means that the number of stable q-ESSs that can be found from the initial 

378 transition is relevant, as this controls how quickly diversity will increase when 

379 a transition does occur in the system. Procedure 2 therefore introduces a high 

380 stochastic variation resulting in a (sometimes sharp, sometimes smooth) di- 
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381 versity increase after an initial (possibly very long) wait. 

382 Procedure 1 bears some resemblance to the colonisation of a new area of land 

383 by many species simultaneously. It results in an initially high diversity as 

384 different q-ESS states form at all points. This decreases quickly to an similar 

385 level found from procedure 2. However, after this time, the two procedures are 

386 equivalent; hence in our analysis we shall consider only initial random seeding, 

387 i.e. procedure 1, in order to standarize the initial diversity level. We then allow 

388 the system to evolve for a long time (40000 generations) before observation to 

389 allow an ecology to form. 




Time (Generations) 

Fig. 2. Wildtype diversity against time for 2 initial systems consisting of the same stable q-ESS at all 
lattice points (initial condition type 2). Diversity remains constant for around 20 generations, after which 
an increase is seen. In one run (solid line), the increase occurs very rapidly but in the other (dashed line) 
the increase is more gradual yet reaches higher levels. Once a stable level (on ecological timescales) of high 
diversity is found, the evolutionary dynamics occurs in the same manner as initial condition type 2, random 
seeding. 

390 The introduction of space has many implications for the model. In the non- 
391 spatial case, there were two timescales: the average lifetime of an individual, 

392 and the average lifetime of a q-ESS, which increased slowly with time. In the 

393 spatial case, we have a third timescale: the time taken for information of a 

394 transition to be transmitted to the other side of the system. As this occurs 

395 only through transitions at all intermediate lattice points, this can be very 

396 long, much longer than the simulation time. Another complication is that 

397 average q-ESS lifetime now depends strongly on the state of neighbouring 

398 lattice points, as migrants from different q-ESSs are disruptive but migrants 

399 from similar q-ESSs can actually stabilise a possible transition. Thus, time 

400 averaging is not possible in large systems, and collecting data on the SAD 

401 becomes very difficult. We therefore focus on calculating the SAR: that is, 

402 the relationship between the number of species found in an area and the size 

403 of the area. We distinguish between the two size measures: the scale as the 

404 sub-area measurement of a system with size X. 

405 SARs come in many forms, depending on the measuring system used. Specif- 
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406 ically, quoting [1], there are 3 main properties : "(1) the pattern of quadrats 

407 or areas sampled (nested, contiguous, noncontiguous, or island); (2) whether 

408 successively larger areas are constructed in a spatially explicit fashion or not; 

409 and (3) whether the curve is constructed from single values or mean values" . 

410 We obtain nested, successive, mean value data. Thus, for all scales, measure- 

411 ment squares are contained within a larger scales' measurement square, no 

412 shapes other than square are considered and we are averaging over all possible 

413 measuring squares from a specific scale. [1] and [2] discuss the implications for 

414 this. 

415 Approximate SAR power-laws are often encountered in real systems at 'medium' 

416 scales: that is, for areas that are smaller than the continent/land-mass that 

417 they are found on, but large enough to obtain a reasonable sample. Good 

418 examples are plant species in Surrey, UK, ([3], page 9) or bird species in the 

419 Czech Republic [10]. When looking at other scales different SARs can be ob- 

420 tained; the distinction between scales is one that varies with environment and 

421 habitat types, and many functional forms of SAR can be found somewhere. 

422 A general rule (p277 of [3]) is that inter-provincial relations follow power-law 

423 SARs with exponent larger than intra-provincially; islands inside a province 

424 will also have a larger exponent than the whole province itself (thus having 

425 smaller diversities). A single run in our model corresponds to a single isolated 

426 province as it is spatially homogenous and self-contained. 

427 A specific instance of our model will not have any real world equivalent, as 

428 we have selected genotype space interactions and our initial position in it 

429 randomly. However, averages over our model should correspond to (large and 

430 thus self averaging) real systems for which our assumptions are approximately 

431 valid, as we are effectively averaging over the possible realisations of genotype 

432 space. Any real world system that does not conform to this average will be 

433 affected by an effect not modelled here - for example, the geography or resource 

434 distribution may be an important factor. 

435 Real systems have z-values between 0.15 and 0.4 depending on the details of 

436 the system [3]. Figure 3 illustrates real SAR data from Hertfordshire plants 

437 and shows a sample simulation SAR. Both describe a power-law as are they 

438 are linear in log-log space, logS 1 = zlogA + logo;, hence the slope of this line 

439 (the z-value) is the major controlling factor in how quickly diversity grows 

440 with area. For example purposes, we have chosen the area of a lattice-point 

441 arbitrarily as 0.4ha. However the true size of a lattice-point in our model is 

442 not well defined as the TaNa model implicitly assumes all species are of equal 

443 spatial extension. Hence we are now concerned only with the scaling relation: 

444 the form of the SAR being close to a power-law and the value of the exponent 

445 in that power-law. 

446 As each run is a separate instance with its own evolutionary history, the diver- 



13 



CO 10 



« Real Data 

— Real data Regression Line (z=0. 19) 
o Simulated data 





Area (ha for data, arbitrary for simulation) 

fa) 



Area (Lattice Points) 

(b) 



Fig. 3. (a) SAR Data for Hertfordshire plants taken from [3] (Fig 2.2) plotted with simulated data, 
assuming 1 lattice-point is a 0.4ha plot (pmove = 0.025) evolved for 40000 generations, (b) Simulated, 
evolved SAR plotted for varying p move from 0.001 to 0.009 (in steps of 0.002); the shape and start point 
remains the same, with only the exponent changing. 

447 sity and z-value variation between runs is high unless the size is much larger 

448 than the species range; however, the power-law rule holds for all instances. 

449 The simulated data in Figure 3 has a slightly reduced tail from the expected 

450 power-law values, due to the finite area of the simulation. By holding a fixed 

451 system size (X = 10 is chosen as be the maximum we can simulate with 

452 sufficient averaging ability) and varying p move (Figure 4 (a)) we can understand 

453 these cutoff's more fully. 
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Fig. 4. (a) z-value calculated from the wildtype diversity evaluated between 40000 and 50000 generations, 
showing individual z-values from runs (on a 10x10 lattice). Note the two distinct regions - p m ove < 0.01 
where species do not spread large enough distances for finite size effects to matter, and p m ove > 0.01 where 
in some runs, species can span the entire system, (b) log-log plot of diversity as a function of time for a 
20x20 system with p m ove = 0.005. 

454 Figure 4(a) shows the individual values of z for varying values of p m ove together 

455 with the average. The values are distributed about some mean, which decreases 

456 approximately linearly with increasing p move for p move < 0.01. however, above 

457 Pmove = 0.01 we observe that some of the runs give a near-zero z- value, i.e. 

458 a constant SAR curve, meaning that species are spanning the system. The 

459 correlation length of the system has reached the system size and boundary 

460 affects will irrevocably effect the results. With increasing p move the average 
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461 patch size of each q-ESS increases, and thus the probability of finding a patch 

462 the size of the system increases. In non-evolutionary models, one can avoid this 

463 problem by considering migration from a 'pool' of constant species makeup 

464 [27] but in evolving systems the pool must be modelled explicitly. 

465 Figure 4(b) shows the time dependence of diversity. Although new species 

466 are produced at all times, and new q-ESS states can be formed, they do not 

467 seem to do so at a rate that matches diversity loss. The time taken to reach a 

468 single q-ESS state diverges with area, taking of the order 10 12 generations for 

469 a single q-ESS to be reached for a 20x20 system, or 10 9 generations for a 10x10 

470 system. As diversity can increase drastically at any time if a single species can 

471 destabilise the dominant q-ESS, it is unlikely this would not continue forever. 

472 Instead, we would effectively be restarting the system with a procedure 2 initial 

473 condition; however, the stability of this highly evolved q-ESS is much higher 

474 than a random q-ESS taken from initial conditions, and so the time taken 

475 to see a restarted system may be very long (as q-ESS lengths are power-law 

476 distributed, this time has mean infinity - however, it does occur eventually, as 

477 there is no truly stable state in this model). 

478 In the Spatial TaNa model, illustrated in Figure 5, the spatial distribution of 

479 species is confined to a contiguous patch. Non-contiguous patches seem to be 

480 rare as patches are more easily invaded at patch corners due to the positive 

481 self-reinforcement of a q-ESS type in the centre. Species will generally coexist 

482 with a specific set of other species, forming fairly distinct q-ESS states of 3-8 

483 species (shaded regions). However, there are many cases where the majority of 

484 q-ESS members remain constant but one species is swapped out for another. 

485 Thus in some cases there is a smooth transition spatially between one q-ESS 

486 type and a completely different q-ESS type, with many transients along the 

487 way containing subsets of each (e.g. dense forest fading to woodlands then to 

488 grassland). In other cases the coexistence is more essential and there will be 

489 a distinct line between one species set and another. 

490 In toroid geometry, any observations of greater than half the total size are af- 

491 fected unaccountably by the periodic boundary so we restrict conclusions from 

492 scales less than X/2, which do appear to be truly power law related (tested for 

493 up to X = 20). Unfortunately, this size restriction does not permit the testing 

494 of self-similarity by any other means than the power-law relation, and we can- 

495 not tell if non-contiguous patches patches might form in larger simulations. It 

496 is possible that species distribution is truly self-similar in our model, whether 

497 the patches are or not, as species may survive in several different patches. 

498 We can also consider this system in the absence of mutation, so considering a 

499 'population dynamics' version of the model. Here, initial conditions are very 

500 important as no new species can ever be added. The quasi-stability observed 

501 previously will also change nature as the only possible disturbance is migrant 
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Fig. 5. Spatial distribution of species on a small (5x5) periodic lattice after 50000 generations, with 
background shading for each point representing the basic q-ESS members and symbols representing all 
genotypes that do not completely fit into a q-ESS category. Some of these genotypes are active in more 
than one q-ESS state (e.g. black circle) and others operate in subsets of a specific q-ESS state (e.g. grey 
triangle). All species are located in contiguous lattice-points, and it is possible for some patches to span the 
entire area. 

502 species. If we for the moment consider a single lattice site with randomly 

503 chosen species, the behaviour is similar to the usual case with mutation in 

504 that the number of species condenses down to a small number which are 

505 mutually stable. As there can be no invasion, the only pressure is accidental 

506 death. This occurs with very low probability for moderate population numbers 

507 as the form of p Q ff ensures that there is a restoring pressure to the equilibrium. 

508 The system will always find a steady state (which, rarely, may have only one 

509 species in if the species that survived the low population stage happen to all 

510 have non-mutualistic interactions). 

511 However, on a spatial lattice things are different. If we choose to evolve a q- 

512 ESS to copy to all points then clearly the system will contain only this q-ESS 

513 forever, as there is no source of change. If we start the system with random 

514 individuals, however, then the initial states found in each lattice point will be 

515 very different and so migrants may have significant impact. In this case, we 

516 see a relaxation in diversity of similar form (power law) to the mutation case. 

517 However, the rate of decay (the exponent for the decrease of diversity with 
sis time) is smaller compared to the evolving case. A species area relation of the 

519 same form as in the evolving case is still seen, complete with slight S shape 

520 form. If we start with an evolved system with a reasonable SAR, and then 

521 turn off evolution, we see that the decay with time of the diversity decreases 

522 drastically, as the system almost 'freezes' (Figure 6). The SAR form will not 

523 change drastically, but the exponent will continue to decrease very slowly as 

524 the number of species, and the number of distinct q-ESS decreases. 

525 This behaviour shows that it is population dynamics that give the SAR power 

526 law form, and that our formalism does not permit mutations to spread through 

527 the system with sufficient speed to offset extinctions. Instead, evolution per- 

528 mits the generation of 'better' q-ESS that can spread through the system 
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Fig. 6. Time dependence of diversity: for the first 5000 generations, mutations are permitted 
(Pmut = 0.001), and are then stopped (averaged over 20 runs). The system decay rate decreases markedly, 
but still follows a power law. 

529 more quickly, accelerating the rate of species loss. However, evolution is re- 

530 quired to produce diversity in the first place, and allows it to spread very 

531 quickly throughout the system as seen in Figure 2. In our model, environ- 

532 mental factors (changing in space and/or evolutionary time) are necessary for 

533 preventing the collapse of the SAR once it is formed. 



534 4 Discussion 



535 Our SAR results bear striking similarity with those of a neutral 'voting' model 

536 of Durrett and Levin [5] . The form of the SAR in both is almost power-law, 

537 with a slight s-shape produced by boundary effects. They find that the z- value 

538 decreases with decreasing speciation rate (which is equivalent to immigration 

539 rate, if new species are introduced from another land mass, for example). In 

540 our model with interactions and explicit genotype space, we find that z-value 

541 decreases with increasing migration rate inside the system. Mutation occurs 

542 at constant speed, so increasing migration rate, e.g. Figure 4(a), decreases 

543 the relative spread of a new species, instead causing transitions to an already 

544 existing q-ESS and so reinforcing currently existing species. 

545 Essentially, internal migration rate reduces the relative effect of mutations, 

546 and so produces the inverse effect of the immigration rate of new species from 

547 outside the system (which is equivalent to mutation in a point-mutation rep- 

548 resentation without consideration of genetics). High mobility (i.e. migration 

549 and immigration rates) for a family of species mean better mixing and so less 

550 chance for spatial segregation of species within a single family - the standard 

551 explanation for why birds generally have lower z-values than land species. 

552 Conversely, e.g. on islands, it allows species from elsewhere to arrive, so possi- 

553 bly increasing diversity (as argued in [5]). Which effect dominates will depend 

554 on the geography in question - i.e. the size of the local groups of individuals, 
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555 and the separation between them. A more detailed model is required to probe 

556 this more fully. 

557 Magurran and Henderson [28], noted that permanent fish species have log- 

558 normal SAD whilst transient species have a log-series distribution. Our local 

559 q-ESS has the same distribution, with a log-normal like distribution for the 

560 wildtypes and a log-series like for mutants and migrants. For low mutation 

561 rates and high migration rates, clearly migrants will outnumber local mutants 

562 and we will observe the exact same distribution near the q-ESS patch borders. 

563 Here, the distinction between the two types is of fitness - the wildtypes with 

564 a log-normal like SAD are all equally fit in that they have a reproduction rate 

565 exactly balancing the death rate; the migrants with a log-series like SAD are 

566 all less fit and rely on repopulation from an external pool. 

567 The Tangled Nature model on a spatial lattice reproduces many of the ob- 

568 served features in real systems without making any a-priori assumptions about 

569 the existence of species. Instead, species and their spatial distributions are al- 

570 lowed to form naturally by co-evolution from simple rules applied only to 

571 individuals. Unfortunately, the model is currently too computer intensive to 

572 allow simulation of the very large scales (and higher migration rates) expected 

573 in real systems. However, a near power law is clearly produced as a simple 

574 result of species forming patches of many sizes, themselves the product of 

575 diffusive dispersion with reproduction and mutation when local interaction is 

576 permitted. Mutation is necessary to give 'raw material' for new species to be 

577 formed. 

578 Co-evolutionary forces are sufficient to allow (co-evolutionary) habitat differ- 

579 entiation (as shown in the co-habitation of competing E.coli strains in [29]), 

580 and the number of different habitats increases with area as a power-law. Thus 

581 power-law SARs are observed, as the number of habitats can drive the diver- 

582 sity increase with area [3], and these persist over long timescales and in the 

583 absence of geographical differences. The evolutionary history therefore relates 

584 to the production, and z-value, of power-law like SARs and may be important 

585 in many cases [3]. 

586 The habitat differentiation produced by co-evolution allows species to be lo- 

587 cally equivalent whilst interacting strongly, and maintains differences in off- 

588 spring probabilities when removed from its favoured habitat. Thus we find 

589 equivalence whenever individuals have had time to adapt to the homoge- 

590 neous killing probability, which corresponds to a situation where individuals 

591 die mainly due to some more our less species independent stochastic killing 

592 mechanism. An example of such a system might be 'climax' stage of forest 

593 succession [30] [31], where species makeup is approximately constant (over a 

594 sufficiently large area and time average) and the ratio of births to deaths are 

595 close to unity for all species. Species measured in the field that were found 
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596 to be non-equivalent [7] may be considered in the context of Tangled Nature 

597 to be transitionary, or may simply be out of the habitat they were originally 

598 adapted to - the equivalence predicted in our system is very local, but can be 

599 formed over distances by the correct migration composition of species. 



eoo Individuals from species not found locally are generally poorly adapted to 

eol the local environment and go quickly extinct. Rarely, however, species with 

602 p f f > can invade and their increased chance of survival over the general 

603 population allows the species to flourish initially - providing a method for fast 

604 speciation from an initial mutant. In addition, during transitions, intermediate 

605 genotypes are successful which may be replaced by other genotypes before a 

606 q-ESS is established, overcoming the 'fitness barrier' to distant genotypes, 

607 with all intermediates occupying fitness maxima. Thus, speciation can occur 
60s quickly, and to species distantly related. This contrasts the 'fitness landscape' 

609 viewpoint (For a review, see e.g. [32]), in which speciation requires passing 

610 through a fitness minima. It also solves a problem seen in neutral theories, 
en which require external pressure such as allopatric speciation (i.e. isolating a 
612 whole community for mutation by "random fission" [33] [34], instead of using 
en the traditional point mutation used here and in much of the literature) if 
en realistically fast speciation and extinctions are to occur [7]. 



615 We have identified the stability of species, fast extinctions and separation in 

616 genotype space as the main differences between our interacting model and 

617 neutral models. The wildtypes in our system are locally equivalent, and it 

618 is the patches of these wildtypes that are producing the power-law SARs 

619 observed. Wildtypes are thus equivalent most of the time but not when found 

620 outside their own habitat, where they suffer a reproductive disadvantage. This 

621 is consistent with the non-neutrality observed in nature and may explain why 

622 neutral dynamics do so well at predicting SARs and SADs. The non-neutrality 

623 is only important during transitions (which, in the spatial model are usually 

624 local events), but the number and distribution of species does not change, only 

625 the specific type of species. These effects cannot be observed in instantaneous 

626 measures, or in time averages. 



627 The spatial Tangled Nature model provides a simple general framework con- 

628 taining the basic properties of diffusive dispersion, reproduction and mutation 

629 on the level of individuals, it allows taxonomic structures to emerge and pro- 

630 duces a large number of observed macroscopic ecological phenomenon - species 

631 abundance, long-lived species, fast extinctions, power-law lifetimes, intermit- 

632 tent dynamics, and, as demonstrated in the present paper, species-area rela- 

633 tions. 
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